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^ Abstract 

cn 

A mutualism is an interaction where the involved species benefit from each other. 
We study a two-dimensional hexagonal three-state cellular automaton model of a 
two-species mutualistic system. The simple model is characterized by four param- 
eters of propagation and survival dependencies between the species. We map the 
parametric set onto the basic types of space-time structures emerged in the mutual- 
istic population dynamic. The structures discovered include propagating quasi-one 
dimensional patterns, very slowly growing clusters, still and oscillatory stationary 
T— I localizations. Although we hardly find such idealized patterns in nature, due to 

1^ increased complexity of interaction phenomena, we recognize our findings as basic 

QQ spatial patterns of mutualistic systems, which can be used as baseline to build up 

<0 more complex models. 
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1 Introduction 



The term symbiosis was coined in 1879 by Heinrich Anton de Bary, a Ger- 
man mycologist, who defined it as the living together of unlike organisms [5]. 
In biology, symbiosis includes all kinds of close relationships between species, 
ranging from parasitic interactions to mutualisms. In two-species interactions, 
mutualism is an interaction of species of organisms that benefits both [36] . 
Most classical symbioses were initially described as two-species interactions, 
but modern research shows that many of these cases comprise more interact- 
ing species and a higher level of complexity than previously thought In 
complex symbiotic systems, the distinction between beneficial and detrimental 
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interactions among the partners is sometimes not clear or can be indirect. The 
outcomes of symbiotic associations may depend also on physiological circum- 
stances, even when only two interacting partners are considered. Mutualism 
is the most intriguing and still not well understood type of inter-species inter- 
actions pi37] . Even simplest models of mutualism exhibits higher behavioral 
complexity than predatoring, parasitism, amensalism and comensalism [3] . In 
the present paper we decided to leave complex interactions and multispecies 
associations beyond the scope, focus on bipartite systems for the sake of sim- 
plicity. This simplification allows us to study idealized localization patterns 
to understand basic laws of spatial structuring without nuisance of real-life 
complexity. 

Mathematical and computer studies of mutualism are growing extensively last 
years. Most analytical models are based on modifications of Lotka-Volterra 
model with positive inter-species interactions [MEj, stabilized with feedbacks 
of limited resources j28j, tuned by diffusiveness and transport effects [17], and 
other relations between interacting species (35|T5] . Other mathematical models 
are based on limit per capita growth and feedback delays 

Almost no results are obtained in spatial evolution of mutualistic popula- 
tions. Spatially extended prey-predator systems produce characteristic wave 
patterns, what are the patterns emerging in mutualistic systems? So far only 
existence of stationary localized domains, or patches, of species is known. 
Their existence is demonstrated by two different techniques: lattice Lotka- 
Volterra model [39] and react ion- diffusion model of population dynamics [3l] . 
Are patches the only patterns which could be observed in mutualistic popu- 
lations? We aim to answer the question in present paper. 

In computer simulations we employ cellular automata: two-dimensional arrays 
of final-state machines which update their states simultaneously and depend- 
ing on states of their immediate neighbours. Cellular automata have been 
used to simulate population dynamics for a long time. They become popular 
as population models from Dewdney articles on prey-predator systems [19], 
and further supported by high-profile research on lattice-gas automata [TU] 
and automata models of host-parasite interaction [27]. For (dis) advantages 
of cellular automata models see overviews [22] and [I6]. The automata mod- 
els are now uncontested models for studying pattern formation in popula- 
tion dynamics [18] , pattern-oriented ecological modeling [26] and spatial ecol- 
ogy [6.40,381^. 

Cellular automata 'substrates' are proved to be successful in imitating and 
simulating propagation of species [TT] , developments of plant populations [1] , 
prediction of epidemics dynamics in spatially heterogeneous environments [21] , 
stochastic species invasion [1^112^ . predation chains [3J and competition in 
complex landscapes [13], and prey-predator systems [H]. Apart of particular 
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case of interacting lattices [9] we are unaware of any published results con- 
certing pure cellular automata models of mutualistic systems. A classical, in a 
sense of Ulam and von Neuman spirit, cellular automata model will be offered 
in present paper. 

In the paper we utilise our ideas on automata-based modelling of population 
dynamics [1], designs of cell-state rules covering all types of inter-species in- 
teractions [2], and our scoping experiments on phenomenology of automata 
model of two-species populations [3]. 

The paper is structured as follows. In Sect. [2] we introduce a cell-state tran- 
sitions rules, imitating mutualism. Phenomenology of the models for various 
parameters of propagation and survival dependencies are outlined in Sect. |3j 
Section |4] discusses basic types of localized patterns discovered and their para- 
metric mappings. 



2 Automaton model of mutualism 

We study hexagonal lattice of finite-state machines, or cells, which take finite 
number of states and update their states in discrete time depending on states 
of their closest neighbours. Every cell x has six neighbours, which determine 
the cell x's neighbourhood u{x) and takes three states: 0, 1 and 2. States 1 
and 2 represent species '1' and species '2'. State represent an 'empty space', 
or a substrate. 

Two processes must be simulated: propagation of species and survival of 
species. A cell of the hexagonal cellular automaton can be occupied exclu- 
sively by 'empty space', state 0, or by one of the species, states 1 or 2. A 
cell X in 'empty state' at time t (x* = 0) becomes occupied by species 1 
at time t + 1 if it has more neighbors in state 1 than in state 2 (cri > a2, 
(Tj = \{y G u{x) : = i}\, i G {1,2}) but enough neighbors in state 2 to 
support species 1 depending on them ((T2 > 6^01), or it has equal amount of 
both species (ai = (72) but there are more species 2 to support species 1 than 
species 1 to support species 2 (a2 — ^^oi > cri — 6*02)- Propagation of species 2 
can be discussed similarly. 

A cell in state 1 (2) remains in state 1 (2) if number of its neighbors in state 
2 (1) exceeds specified threshold, a2 > 6*11 (cti > 6^22) ■ There are no cell-state 
transitions 1^2 and 2^1. See summary of the cell-state transitions in 
Fig.[Tj 

The thresholds 6'oi and 6q2 are propagation, or sustainability parameters. 
The thresholds On and 6*22 are survivability parameters. In general, 1 < 



3 










1 


2 





(71 + 


0-2=0 




cy {^c AD) 


1 


0-2 


<^11 


(72 > 6'ii 


never 


2 




< ^22 


never 


(Ti > 6122 



Fig. 1. Cell-state transition table imitating mutualistic in- 
teractions. The conditions A, B, C and D are as follows: 
A = {ai > 0-2) A (cJ2 > 6*01) V (o-i + (T2 > 0) A {a2 - Oqi > ai - 602), 
C = (cJ2 > fJi) A (cti > 6*02) V (cJi + a2 > 0) A (cJ2 - 6^01 < cxi - 6^02), 
B = (random(l) < 0.5), D = (randoni(l) > 0.5), and random generates a random 
real number between and 1. 

^oi! ^02, 6'22 < 6, however no activity persists on the lattice for values higher 
than 3, so in the paper we consider only 1 < 6*01, 6*02, 6*11, 6*22 < 3. 

Sometimes we are addressing a cell-state rule as rule R{9oi, 6*02, 6*11, 6*22)- 



3 Phenomenology 



We have analysed development of automata from random initial configurations 
for all values of four thresholds 1 < 6'oi, ^02, ^11, 6*22 < 3. We found that when 
one of the species 1 (2) does not depend on another species to propagate 
6*02 = (6*01 = 0) then this species 1 (2) propagates unlimitedly while other 
species remains confined to the zone of initial inoculation. See comprehensive 
list of configurations in Appendix and few selected samples in Fig. [2j 

Dynamics of unlimited propagation depends on parameters 601 and ^02- Con- 
sider for example, 6'oi = and 6*02 < 4. Species 1 propagates on the lattice 
because they do not need species 2 to invade new space. If 9u = then lattice 
becomes filled with 'solid' pattern of species 1, once a site occupied by 1 it 
always remains in the state 1. For On > we observe structures similar to 
target waves because species 1 can not remain at the same site for more then 
one step of discrete time t, it dies, or vacates the site, at time t + 1, however 
the site becomes occupied again by species 1 at time step t + 

For the parameters ^02 > and ^01 > localizations — compact groups of 
non-quiescent (not equal 0) states — are observed. The number of localizations 
in configurations, recorded after transient period, decreases with increase of 
^02 + ^01- Principle localizations are discussed in next section. 
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Fig. 2. Samples of configurations parameterized by values of ^oi; rows, and 602, 
columns. Each entry [a, /3] in the table represents a sample configuration for pa- 
rameters ^01 = 9q2 = P, Oil = 1 and 622 = 2. Lattice is 200 x 200. Initially 
cells not further then 50 cells from the lattice centre are assigned states 1 or 2 with 
probability 0.1 each. 

4 Localizations in mutualistic systems 



In exhaustive experiments we discovered four principle types of localizations: 
stationary still localizations, stationary oscillating localizations, propagating 
quasi-one-dimensional localizations, or worms, and sub-linearly growing do- 
mains. 'Stationary' means that the pattern does not travel along the lattice, 
'still' means that the pattern does not change its structure during automaton 
development. 

Stationary still localization is a compact groups of cells in non-quiescent (non- 
'0') states. Diameter of such group usually does not exceed 1-3 cells. Typically 
one species occupies central position in the localization and it is surrounded by 
sites occupied by another species. Thus in Fig. [3] we see sites of species 1 (solid 
discs) surrounded by several sites with species 2 (circles). Species 1 survive 
because they have more than two neighbouring sites with species 2 (6*02 > 2) 
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Fig. 3. Stationary still localizations in rule i?(2221). State '1' is shown by solid disc, 
state '2' by circle, state '0' by dot. 



Fig. 4. Oscillators in rule i?(2323). 
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Fig. 5. Oscillators in rule i?(1223). 

while less dependent species 2 are supported by being in relation with just 
one site occupied by species 1 (^oi > !)• Several elementary localizations can 
be linked together to form an extended cluster of localizations. 
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Stationary oscillator is compact group of cells in non-quiescent states, which 
does not move along the lattice unlimitedly. The group can be translated cycli- 
cally around some fixed point. All oscillators observed in space-time dynamics 
of mutualistic populations are can be classified as fiip-fiops, which rotate on a 
fixed degree around stationary centre, e.g. oscillators in Fig. |4] and Fig. |5} and 
breathing oscillators, which has stationary core, e.g. sites in state 1 in Fig. |6} 
and switching/breathing halo, e.g. sites in state 2 in Fig. [6] 

Some rules support only still localizations, e.g. rule i?(2221) (Fig. [3]), some only 
oscillators, e.g. rule i?(2323) (Fig. [i]), others support both still and oscillating 
localizations, e.g. rule i?(1113) (Fig. [6]). 

Worms are quasi-one-dimensional propagating localizations, which usually 
consists of two parallel chains of non-quiescent states, one chain is formed 
by sites with species 1, another chain by sites with species 2. All worms ob- 
served in our experiments have only two growing tips. For example, see in 
Fig. [T] snapshots of a single worm growing in rule i?(llll) automaton. The 
worm has two growing tips which are two sites: one occupied by state 1, an- 
other by state 2. For neighbouring quiescent cells numbers of non-quiescent 
states 1 and 2 is equal, so the states 1 and 2 'diffuse' randomly. Species 1 
and 2 can not survive without each other, therefore when one of the species 
advances, it must then 'wait' for a neighbouring empty site to be occupied 
by another species. Due to randomization involved the growing tips do not 
propagate directly but rather implement a kind of random walk on the lattice 
(Figgi). 

When initial density of species 1 and 2 is high enough, many worms are born 
(FigjsjD), they propagate and collide with each other. When a growing tip 
collides to a body of a worm it stops propagating. Thus forming a porous 
stationary core with free growing tips propagating centrifugally (FigjsjD). 

For some rules, e.g. rule i?(3300) (Fig.|9]), worms are short-living. In these rules 
worms stop propagating when their growing tips are blocked by the following 
conditions. Every resting cell x, = 0, closest to non-resting cells of the 
worm's tip has {a\ > al) A {al < 9qi), i.e. species 1 do not propagate, and also 
(o"2 > al) A (a\ < 602), i.e. species 2 do not propagate as well. 
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Fig. 6. Still localizations and oscillators in rule -R(1113). (a)-(c) three consecutive 
configurations exhibiting still localizations: sites with species 2 surrounded by sites 
with species 1, and stationary oscillators: sites with species surrounded by sites with 
species 2. (d) zoomed out configuration. 
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(a) f = so (f) * = SI (g) t = 52 th) t = 

Fig. 7. Propagation of long-life worms in rule i?(llll). 
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(a) (b) 

Fig. 8. Long-life worms in rule i?(llll). (a) a single worm grown for 197 steps 
of discrete time, (b) Zoomed out configuration of rule i?(llll) automaton, which 
starts its initial development in a random configuration, where all cells closer then 
150 sites to the centre of the lattice are assigned states 1 or 2 with probability 0.5. 
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Fig. 9. Propagation of short-life worms in rule -R(3300). 
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Fig. 10. A snapshot of a sub-linearly growing stripy domains, rule i2(2222), taken 
at llTlth step of automaton development from initial random configuration. 



The last class of localizations are very slowly growing stripy domains (Fig. [Io|), 
found in space-time dynamics of cellular automata governed by rules i?(2222) 
and i?(3322). They are not proper localizations, because their size is con- 
stantly increasing. However it does increase very slow: from computational 
experiments we estimated that a radius of a circle around a stripy domain 
grows as 0.05 ■ t, where t is a number of time steps from the beginning of the 
pattern development. That is the growth rate is almost two orders smaller then 
typical propagation rate of usual cellular automaton patterns (e.g. excitation 
waves or gliders in Conway's Game of Life). So comparing to typical cellular 
automaton growing patterns the stripy domains behave as localizations. 

The domains are stripy, i.e. consist of intermittent chains of 1 and 2 states, 
because each species needs at least two or three sites occupied by other species 
to propagate and at least two sites occupied by other species to survive. Thus 
having a chain of states 1, the chain will be covered by on both sides by states 
2, which will be covered by states 1 and so on. See mechanics of the pattern 



formation in Fig. 11 Growth of the patterns is slow due to involved of random- 
ization, see conditions B and C in cell-state transition table (Fig. [T]). When 
an empty site accept species 1 or 2 at random, the species may not stay there 
if survivability criteria not met, then the site become empty again and may 
take of the species at random. For example, see how species 2 propagate to the 
North tip of the growing pattern in Fig. [TT^-x. At first empty site is occupied 
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Fig. 11. Mechanics of stripy pattern formation in rule i?(2222). 
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by species 2 (also because there are neighbouring species 1), however species 
1 disappear because their neighbourhoods do not meet survivabihty criteria; 
and the site just occupied by species 1 becomes empty again (Fig. 11 *c). 



5 Conclusions 

In the paper we introduced a minimal cellular automaton model of a two- 
species mutualism, where every cell of a two-dimensional hexagonal lattice 
takes three states: species 1 and 2 and empty space 0, and updates its state 
depending on number of another species in its neighbourhood. The cell-state 
transition rules are threshold-based: a species propagate if number of another 
species exceeds certain threshold (propagation dependence) and the species 
continue occupying the newly acquired site if number of another species exceed 
certain threshold (survivability dependence). We were mainly concerned with 
a particular type of mutualism when both species can not propagate without 
each other but can survive without each other. 

We have discovered four types of localizations — loci of non-empty sites which 
either preserve their size or grow very slowly. They are stationary still localiza- 
tions, stationary oscillators, quasi-one-dimensional propagating domains and 
sub-linearly growing stripy domains. 

We were looking in the natural world for localization which resemble the ones 
presented in this work. However, we were unable to find any of the patterns 
in natural symbiotic systems. Evidence for spatial self-organization is present 
in mussel communities. Small-scale, labyrinthine spatial patterning, which re- 
semble somehow our sub-linearly growing stripy domains, are present in beds 
of the blue mussel Mytilus edulis [BOJ. However, these patterns are not the 
results of symbiotic interactions. Previous work showed that an empirically 
derived cellular automaton approach could predict macroscale patterns |12]. 
The rather complex full cellular automaton model comprised 15 ecological 
states which also include states for co-occurring other species. The spatial 
structures formed by this model are patchy clusters of mussels with separat- 
ing gaps, with their size distributions being dependent on the inclusion of 
interactions with other species. 

Why do we not find the idealized patterns in nature ? Even the simplest 
cases of inter-species interactions are apparently more complex. Many species 
produce diffusible compounds which can have signalling or nutritional func- 
tion, blurring the local spatial consequences of interactions among the species, 
which and can change the behaviour of species, such as exemplified by quorum- 
sensing. Colonization of natural substrates is usually not restricted to two di- 
mensions. In early or late stages bacterial multispecies biofilms grow in the 
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Propagating worms: J?(1100) , 
K(llOl), i2(llll); Very slowly- 
growing worms / domains: Ji(1122) ; 
Still and oscillating domains for 
other values of parameters. 



All clusters arc non-propagating. Ma- 
jority are stationary still clusters. 
Only for /?(1222), H(1223), fi(1233) 
oscillating clusters. 

Worms propagate but freeze quickly: 
H(2200), H(2201), J2(2202); For 
H(2202) life-time of worms is really 
short and no worms appear for 
il(2203); for ii(2211) worms freeze 
at once and form still clusters. For 
H(2212) and H(2213) still clus- 
ters, J?(2222) slowly growing stripy 
domains. 



Mostly stationary still clusters, os- 
cillating clusters only for J?{1322) , 
i?(1323), il(1333). 



Mostly still stationary clusters, oscil- 
lating only for H(2322), J?(2323), for 
i?(2333) all activity extinguish all dy- 
ing. 



For /?(3300) worms propagate but 
freeze quickly; for i2(3301)-fi(3303) 
number and life-time of warms 
decrease, till no worms at all for 
H(3303); H(3311) long-living worms; 
for i?(3312) to i2(3313) worms 
quickly freeze or die. For J2(3322) 
worm-like large clusters grow. There 
is no sustained activity for J?(3323) 
and J2(3333) no sustained activity. 



Fig. 12. Characteristics of space-time dynamics for Oqi > 0, rows, and 602 > 0, 
columns. Each entry [a, (3] provides a brief description of space-time dynamics 
of a sample configuration for parameters 6qi = a, 6q2 = /?, < On < 3 and 
< 6*22 = 2 < 3. 



third dimension to produce characteristic shapes and phenomena, which are 
often associated with altruistic production of extracellular polysaccharides in 
one or more of the interacting species [H]. In addition, the substrate is often 
not exactly homogeneous, and even minute variations of surfaces can cause re- 
sponses of the emerging community. Finally, natural symbionts can be motile 
to translocate the partner species. However, our simplified model represents 
principal rules of spatial patterning in mutualisms which may serve as starting 
points to buid up more complex models of associations. 



Space-time dynamic of mutualistic populations parameterized by propaga- 



tion dependencies is shown in Fig. [T2j There we see that when propagation 
dependencies increase or when dependence of one species is higher then depen- 
dence of another species the stationary localizations change from propagating 
worms to oscillation localizations to stationary localizations. Most remarkable 
varieties of localizations emerge for lower but above zero propagation depen- 
dencies. The very similar correlation takes place between spatial patterns and 
survivability parameters (see full catalog of patterns in Appendix), i.e. tran- 
sitions from worms to oscillators to still localizations. 
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